## This code include scripts for Appendix Figures
# Figure A1-a and A1-b

rm(list = ls())

library(readr)
library(ggplot2)
library(stargazer)
library(lme4)
library(dplyr)
library(haven)
library(texreg)
library(tidyr)

# variables:
#"onset_ko_flag":Binary flag indicating group-level conflict onset / ko option
# "onset_ko_terr_flag": binary flag indicating group-level territorial conflict onset  / ko option
# "onset_ko_gov_flag: binary flag indicating group-level governmental conflict onset  / ko option
# onset_do_flag": binary flag indicating group-level conflict onset: / do option
# "onset_do_terr_flag": binary flag indicating group-level territorial conflict onset  / do option  
# "onset_do_gov_flag": binary flag indicating group-level governmental conflict onset  / do option  


load("Final_JPR/data/paperdata.RData")
data <- paperdata %>%
       filter(tek_count>0)

source("Final_JPR/Code/coefplot.R")
source("Final_JPR/Code/marginal_effects.R")
source("Final_JPR/Code/setup.R")


## dependent variable
dv <- 'onset_do_flag' 
## controls

ivs_1 <- c('best_tlineq_total', 'tekbest_status_excl')
ivs_2 <- c('worst_tlineq_total', 'tekworst_status_excl')
ivs_3 <- c("median_tlineq_total",'tekmedian_status_excl')
ivs_4 <- c('nearst_tlineq_total','teknearst_status_excl')


f_do_1 <- formula(paste0(dv, ' ~ ', paste(paste(ivs_1, collapse=' + '), sep =' + ')))
f_do_2 <- formula(paste0(dv, ' ~ ', paste(paste(ivs_2, collapse=' + '),  sep =' + ')))
f_do_3 <- formula(paste0(dv, ' ~ ', paste(paste(ivs_3, collapse=' + '),  sep =' + ')))
f_do_4 <- formula(paste0(dv, ' ~ ', paste(paste(ivs_4, collapse=' + '),  sep =' + ')))

ivs <- c('best_tlineq_total','worst_tlineq_total','median_tlineq_total','nearst_tlineq_total')


## using log of transnational inequality
model_do_log <- list()
for (i in 1:length(ivs)){
       model_do_log [[i]] <- glm(eval(parse(text = paste0('f_do_', i))), 
                                 data = data, family = binomial(link = "logit"))
}



p1 <- marginaleffect_logit(ModelResults = model_do_log, n.sim = 1000, 
                           varname = c("best_tlineq_total", "worst_tlineq_total",
                                       "median_tlineq_total","nearst_tlineq_total"),
                           data = data,
                           val1 = 0, val2 =1, clusterid = "gwgroupid") 

p1 +  scale_x_continuous(breaks = 1:4,
                         labels = parse(text =rev(c( 
                                "Nearest~TEK",
                                "Median~TEK",
                                "Worst~TEK",
                                "Best~TEK")))) + ggtitle("Conflict Onset")
ggsave("Final_JPR/Figures/appendix_figures/fig_A1-a.pdf", width = 6.5, height = 4.2, units='in', dpi=600)


####################################### Model 1 across Tables A1-5 ######################################################################
## dependent variable
dv <- 'incidence_flag' 
## controls


ivs_1 <- c('best_tlineq_total', 'tekbest_status_excl')
ivs_2 <- c('worst_tlineq_total', 'tekworst_status_excl')
ivs_3 <- c("median_tlineq_total",'tekmedian_status_excl')
ivs_4 <- c('nearst_tlineq_total','teknearst_status_excl')


f_do_1 <- formula(paste0(dv, ' ~ ', paste(paste(ivs_1, collapse=' + '), sep =' + ')))
f_do_2 <- formula(paste0(dv, ' ~ ', paste(paste(ivs_2, collapse=' + '),  sep =' + ')))
f_do_3 <- formula(paste0(dv, ' ~ ', paste(paste(ivs_3, collapse=' + '),  sep =' + ')))
f_do_4 <- formula(paste0(dv, ' ~ ', paste(paste(ivs_4, collapse=' + '),  sep =' + ')))
#
ivs <- c('best_tlineq_total','worst_tlineq_total','median_tlineq_total','nearst_tlineq_total')


## using log of transnational inequality
model_incidence_log <- list()
for (i in 1:length(ivs)){
       model_incidence_log  [[i]] <- glm(eval(parse(text = paste0('f_do_', i))), 
                                         data = data, family = binomial(link = "logit"))
}

p2 <- marginaleffect_logit(ModelResults = model_incidence_log, n.sim = 1000, 
                           varname = c("best_tlineq_total", "worst_tlineq_total",
                                       "median_tlineq_total","nearst_tlineq_total"),
                           data = data,
                           val1 = 0, val2 =1, clusterid = "gwgroupid") 

p2 + scale_x_continuous(breaks = 1:4,
                        labels = parse(text =rev(c( 
                               "Nearest~TEK",
                               "Median~TEK",
                               "Worst~TEK",
                               "Best~TEK"))))+ ggtitle("Conflict Incidence")
#ggsave("JPR_RR/Tex/Figures/fd_do_incidence_control.pdf", width = 6.5, height = 4.2)
ggsave("Final_JPR/Figures/appendix_figures/fig_A1-b.pdf", width = 6.5, height = 4.2, units='in', dpi=600)




